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ABSTRACT 


This  thesis  explores  the  numerical  methods  and  software  development  for 
optimal  trajectories  of  a  specific  model  of  Helicopter  Unmanned  Aerial  Vehicle 
(UAV)  in  an  obstacle-rich  environment.  This  particular  model  is  adopted  from  the 
UAV  Laboratory  of  the  National  University  of  Singapore  who  built  and  simulated 
flights  for  an  X-Cell  60  small-scale  UAV  Helicopter.  The  code,  which  allowed  the 
team  to  simulate  flights,  is  a  complex  system  of  non-linear  differential 
equations — 15  state  variables  and  four  control  variables — used  to  maneuver  the 
state  trajectories.  This  non-linear  model  is  incorporated  into  a  separate 
optimization  algorithm  code,  which  allows  the  user  to  set  initial  and  final  time 
conditions  together  with  various  constraints,  and,  using  the  same  variable 
scheme,  optimize  a  trajectory.  The  optimal  trajectory  is  defined  by  using  a  cost 
function — the  performance  measure — and  the  system  is  subject  to  a  set  of 
constraints  (such  as  mechanical  limitations  and  physical  three-dimensional 
obstacles).  Simulations  conclude  that  solutions  are  readily  obtained;  however,  it 
is  still  very  difficult  to  derive  trajectories  that  are  truly  optimal,  and  our  work  calls 
for  more  future  research  in  computational  programs  for  optimal  trajectory 
planning.  All  simulations  in  this  thesis  are  modeled  using  the  MATLAB  program. 
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I.  LITERATURE  REVIEW 


A.  BACKGROUND 

This  project  explores  and  simulates  a  specific  algorithm  for  trajectory 
optimization  for  UAVs  in  an  environment  with  fixed  or  unfixed  obstacles  that  may 
be  hard  or  soft  in  nature.  A  hard  obstacle  represents  a  physical  obstruction 
where  contact  will  result  in  damage  to  the  system  and  mission  failure.  A  soft 
obstacle  represents  something  that  should  be  avoided,  but  the  algorithm  user 
can  designate  a  priority  or  weighting  to  the  obstacle  to  portray  its  level  of 
detriment  to  the  UAV.  Since  both  path  and  trajectory  planning  are  popular  topics 
of  research  among  the  science  community — with  a  plethora  of  applications — it  is 
important  to  conduct  a  literature  review,  both  to  gather  ideas  that  may  be  helpful 
to  research  and  simulation,  and  to  ensure  that  our  algorithm  is  unique  in  nature. 

There  are  two  types  of  relevant  problems  in  the  reviewed  literature:  path 
planning  and  trajectory  planning  and  optimization.  Path  planning  is  the  least 
complex  of  the  two,  but  may  provide  some  simple  ideas  that  can  still  be  applied 
to  more  complex  problems.  Trajectory  planning  and  trajectory  optimization  are 
very  similar  as  far  as  variables  and  parameters  are  concerned.  The  main 
difference  is  that  trajectory  optimization  is  focused  on  maximizing  or  minimizing  a 
specific  cost  function  (i.e.,  travel  time,  cost,  fuel  consumption).  This  chapter 
provides  a  brief  description  of  pertinent  articles  researched  in  each  of  the  two 
categories. 

B.  ARTICLE  REVIEWS 

Path  planning  is  the  simplest  form  of  planning,  since  it  generally  involves 
directional  components  of  a  vehicle  or  aircraft,  which  means  that  the  model  will 
consist  of  two  space  variables  for  2-D  and  three  space  variables  for  3-D. 
Occasionally,  an  analysis  of  path  planning  will  include  velocities  which  would 
increase  the  number  of  variables  to  four  and  six,  respectively.  The  best  path  is 

created  by  optimizing  these  variables  over  a  set  of  constraints. 

1 


In  [1],  a  team  analyzes  a  problem  slightly  more  difficult  than  standard  path 
planning,  as  it  describes  a  method  for  finding  the  shortest  path  for  multiple  UAVs 
flying  simultaneously,  rather  than  just  one  aircraft.  This  team  uses  Dubins  paths, 
which  calculate  a  sum  of  circular  arcs  and  their  tangent  lines,  to  determine  the 
best  paths  for  a  swarm  of  UAVs  to  simultaneously  converge  on  a  target  region. 
The  model  uses  three  general  constraints:  (1)  A  maximum  bound  on  the  UAV 
curvature  (due  to  limitations  of  the  aircraft),  (2)  Minimum  separation  distance 
between  UAVs,  and  (3)  Non-intersection  of  paths  at  equal  length  (to  prevent 
collision).  The  model  intends  to  protect  the  aircraft  while  minimizing  distance 
travelled,  which  inherently  achieves  the  goals  of  minimizing  fuel  consumption 
while  increasing  durability  of  the  UAVs  [1], 

Each  UAV  path  is  generated  through  finding  a  common  plane  between  the 
initial  and  final  position  vectors  and  the  final  tangent  vector  (i.e.,  for  each  path, 
the  three  vectors  are  co-planar).  Then,  each  path  length  is  calculated  using  the  a 
specific  set  of  equations.  Each  path  is  compared  to  a  “Reference  Path”  which  is 
the  longest  of  all  possible  UAV  paths.  Each  path  is  then  lengthened  (by  reducing 
curvature)  and  compared  against  all  existing  constraints  until  an  optimal  solution 
is  obtained.  Since  no  aircraft  specifications  pertain  to  this  simple  model,  it  may 
be  used  for  fixed  or  rotary  winged  aircraft  [1], 

The  aspects  of  this  article  are  helpful,  since  they  provide  a  method  for 
coordinating  several  UAVs  for  simultaneous  arrival  on  a  designated  target; 
however,  they  avoid  several  critical  aspects  that  are  better  addressed  in 
trajectory  planning  and  optimization.  First,  this  model  assumes  an  equal  path 
length  and  a  constant  speed  for  all  UAVs.  Second,  the  analysis  assumes  an 
obstacle-free  environment.  Finally,  the  equations  used  for  this  path  planning 
method  do  not  optimize  any  critical  performance  functions  that,  in  turn,  do  not 
allow  any  maximization/minimization  planning  effects  [1], 
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Models  of  trajectory  planning  and/or  optimization  become  more 
complicated  since  they  incorporate  dynamics  involving  measures  of  performance 
of  the  actual  vehicle  or  aircraft,  as  well  as  directional  variables.  While 
parameters  and  constraints  still  exist  in  the  model,  the  number  of  variables  and 
dimension  of  the  differential  equation  system  will  increase  substantially.  In  this 
case,  the  nonlinearity  in  the  dynamics  becomes  the  main  challenge  in  the 
searching  for  optimal  solutions. 

In  [2],  the  authors  study  a  time-variant  model  by  analyzing  the  probability 
of  a  UAV  becoming  disabled  at  a  certain  time  juncture  in  a  dynamic  environment. 
In  other  words,  the  probability  of  a  UAV  becoming  disabled  varies  over  the  time 
variable.  This  particular  model  requires  that  the  UAV  reach  the  objective 
(accomplishing  the  mission)  while  utilizing  the  shortest  path  possible  or  finding 
the  trajectory  that  minimizes  the  probability  of  the  UAV  being  disabled. 
Weighting  of  parameters  by  the  user  determines  which  criterion  is  most  important 
to  the  model.  The  probabilistic  map  changes  constantly  over  time;  therefore,  the 
paths  generated  by  this  strategy  incorporates  functions  of  both  position  and  time. 
Consequently,  variables  will  exist  for  velocity  and  acceleration,  not  just  for  the 
position  vector  (as  is  the  case  in  path  planning).  That  means  that  this  model 
must  also  incorporate  several  constraints  for  UAV  capabilities.  As  a  UAV  flies  in 
an  area  with  multiple  threats,  the  risk  of  the  UAV  becoming  disabled  is 
characterized  by  the  probability  of  the  UAV  becoming  disabled  at  a  certain 
location  and  time.  The  probability  is  modeled  using  Gaussian  Probability  Density 
Functions  [2], 

This  probabilistic  approach  provides  a  much  more  useful  method  for  UAV 
planning,  using  trajectory  planning  rather  than  path  planning.  This  allows  the 
modeler  to  emplace  constraints  based  on  the  capabilities  of  the  UAV  and 
capabilities  of  the  threat  as  well  as  path  direction.  The  model  also  allows  all 
variables  to  change  over  time  which  means  that  obstacles  can  be  fixed,  moving, 
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and/or  with  changing  capabilities  based  on  both  time  and  position.  These 
aspects  of  this  model  are  very  useful  for  the  prospect  of  real-time  trajectory 
planning  [2], 

In  [3],  CDR  Mike  Hurni  analyzes  optimal  control  techniques  for  vehicle 
trajectories  specifically  pertaining  to  minimizing  the  cost  function  for  ground 
vehicle  operations.  His  dissertation  shows,  very  methodically,  how  optimal 
control  techniques  can  work  effectively  to  minimize  cost  (or  another  user  selected 
performance  measure)  while  maximizing  the  flexibility  for  a  ground  vehicle  to 
alter  its  path  in  a  real-time  fashion  in  response  to  obstacles  known  or  not  known 
a  priori.  In  other  words,  his  techniques  can  minimize  cost  while  maximizing 
robustness.  His  techniques  also  provide  a  great  deal  of  flexibility  in  allowing  the 
user  to  alter  weighting  values  to  increase,  decrease,  or  eliminate  the  impact  of 
certain  aspects  of  the  trajectory  on  the  overall  cost  function  value  [3], 

Additionally,  CDR  Hurni’s  work  provides  a  series  of  “requirements 
checks”  on  a  particularly  selected  trajectory  to  test  feasibility,  obstacle 
collision/buffer  violation,  caution  zone  infractions,  and  several  other  constraints 
that  may  deem  the  trajectory  to  cause  mission  failure.  These  checks  are 
established  in  a  way  to  allow  the  user  to  vary  buffer  zones,  the  number  of 
obstacles,  the  state  of  obstacles — fixed  or  unfixed — and  other  key  criteria  that 
are  critical  to  the  flexibility  of  this  particular  optimal  control  system  [3], 

In  the  final  stages  of  the  dissertation,  CDR  Hurni  addresses  the  prospect 
of  multiple  ground  vehicles  operating  in  the  same  space  simultaneously.  His 
algorithm  allows  the  user  to  implement  weights  involving  collision  and  grazing 
prevention  for  a  varying  number  of  vehicles.  These  weights  are  in  addition  to  the 
previously  calculated  weights  for  obstacle  avoidance.  This  addition  allows  an 
even  broader  range  of  use  for  the  algorithm.  CDR  Hurni’s  use  of  optimal  control 
methods  has  resulted  in  an  extremely  flexible  algorithm  that  allows  a  user  to 
input  and  alter  numbers  and  weightings  for  obstacle  and  ground  vehicle  criteria 
and  constraints  while  optimizing  trajectory  to  minimize  or  maximize  a  function 
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selected  by  the  user.  The  only  disadvantage  to  this  project,  as  it  pertains  to  our 
work,  is  that  our  problem  deals  with  a  nonlinear  system  and  much  more  complex 
dynamics  [3], 

In  [4],  a  team  documents  their  construction  of  a  vehicle  to  compete  in  the 
DARPA  Urban  Challenge  by  formulating  and  solving  an  optimal  and  real-time 
trajectory  planning  problem  with  velocity  variables  and  nonlinear  dynamics.  In 
forming  their  optimal  trajectory  formula,  the  team  integrates  boundary  conditions 
at  the  initial  and  final  time  instances,  motion  constraints,  and  collision  avoidance 
criteria  as  the  three  conditions  that  must  be  satisfied  to  enable  a  feasible 
trajectory.  The  fourth  constraint  is  the  minimization  of  the  “performance  index” 
(optimization).  These  constraints  are  designed  with  the  assumption  that 
boundary  conditions  and  motion  constraints  are  generally  given  while 
optimization  and  collision  avoidance  criteria  are  chosen  to  fit  the  designer’s 
needs  for  a  specific  situation  [4], 

The  key  development  for  this  trajectory  planning  problem  is  the  idea  of 
real-time  obstacle  avoidance.  This  means  that  it  is  assumed  that  obstacles  may 
be  fixed  or  moving  and  that  the  positions  of  these  obstacles  are  generally  not 
known  a  priori  which  means  that  trajectories  must  be  re-planned  multiple  times 
throughout  movement  between  the  initial  and  final  positions.  Three  key  features 
enable  this  method  to  satisfy  the  requirement:  (1)  All  paths  satisfying  boundary 
conditions  and  the  vehicle’s  kinematic  constraints  are  parameterized  in  terms  of 
polynomials  of  sufficient  order,  (2)  A  collision-free  criterion  is  developed  and 
imposed  for  avoiding  both  “hard”  and  “soft”  obstacles  that  are  detected  along  the 
path  (in  real  time),  and  (3)  A  performance  index  is  introduced  to  find  the  best 
path  among  the  collision-free  paths  meeting  all  necessary  criteria  [4],  Finally,  the 
performance  index  is  chosen  so  that  paths  equivalent  to  the  shortest  path  can  be 
solved  analytically  while  meeting  all  criteria — based  on  both  given  constraints 
and  those  imposed  by  the  designer.  The  bottom  line  is  that,  for  a  successfully 
optimal  trajectory  that  completes  the  mission  in  a  real-time  changing  dynamic 
environment,  the  trajectory  must  be  modeled  in  a  piecewise  fashion  as  the 
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vehicle  moves  from  its  initial  to  final  position.  The  team’s  conclusion  is  that  this  is 
most  effectively  accomplished  through  use  of  a  parameterized  fourth-order 
polynomial.  By  obtaining  the  solution  to  an  adjustable  parameter,  it  is  possible  to 
generate  a  real-time  solution  [4], 

The  method  discussed  in  this  article  provides  the  same  advantages  as 
were  discussed  in  CDR  Hurni’s  dissertation  with  the  exception  of  not  including  a 
model  with  multiple  vehicles  to  be  coordinately  controlled.  Also,  the  solution  for 
this  model  was  derived  analytically,  and  our  complex  nonlinear  model  will  most 
certainly  require  a  numerical  model  to  ensure  solvability. 

C.  CONCLUSION 

The  professional  works  uncovered  during  the  literature  review  phase  of 
this  thesis  have  shown  that  trajectory  optimization  is  a  complex  problem  which 
complicates  drastically  when  non-linear  dynamics  and  constraints  are  involved. 
The  papers  reviewed  demonstrate  several  methods  for  analyzing  and  solving 
path  and  trajectory  planning  problems.  They  all  have  strengths  and  weaknesses 
that  are  pertinent  to  our  future  research.  The  key  point  is  that  none  of  them  truly 
solve  our  exact  situation  that  requires  optimality  of  the  performance  measure,  so 
we  will  surely  incorporate  some  of  the  ideas  from  these  works.  Significant 
original  input  will  also  be  required  to  fulfill  our  objective:  to  derive  and  test  an 
algorithm  to  numerically  solve  an  optimal  trajectory  for  a  UAV  in  a  dynamic 
environment  with  and  without  obstacles. 
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II.  PROBLEM  FORMULATION 


A.  INTRODUCTION 

At  the  end  of  Chapter  I,  it  is  clearly  stated  that  the  objective  is  to  produce 
an  optimal  UAV  trajectory  in  a  dynamic  environment  with  fixed  and/or  moving 
obstacles.  This  objective  requires  the  derivation  of  a  system  of  equations  and 
functions  that  consider  three  separate  elements:  1)  The  dynamic  equations  that 
determine  the  trajectory  of  the  UAV,  2)  The  trajectory  constraints,  and  3)  The 
performance  measure  for  optimization.  The  dynamic  equations  involve  a  system 
of  both  state  and  control  vectors.  The  state  vectors  are  represented  by  the 
positions,  attitudes,  and  velocities  of  the  UAV — all  with  respect  to  time.  The 
control  vectors  are  determined  by  the  system  inputs  controlled  by  the  algorithm’s 
user.  Trajectory  constraints,  in  this  case,  alter  the  trajectory  to  avoid  all 
obstacles  in  a  projected  path.  The  optimizing  performance  measure  outlines  a 
function  to  minimize  cost,  time,  or  another  parameter  of  our  choosing  by  solving 
for  the  optimal  control  vector  in  a  given  trajectory.  The  remainder  of  this  section 
will  use  the  three  elements  explained  above  to  develop  a  logical  and  solvable 
problem  that  allows  us  to  fulfill  the  stated  objective. 

B.  TRAJECTORY  EQUATIONS 

The  trajectory  equations  for  this  particular  system  will  be  structured  in  the 
form  of  a  comprehensive  non-linear  differential  equation  model  based  on  a 
Helicopter  UAV  built  and  modeled  by  a  research  team  from  the  National 
University  of  Singapore.  The  model  will  consist  of  four  key  components  that  will 
derive  a  total  of  15  state  variables  and,  consequently,  15  differential  equations: 
1)  Kinematics,  2)  Six  degree-of-freedom  (DOF)  rigid-body  dynamics,  3)  Main 
rotor  flapping  dynamics,  and  4)  Yaw  rate  gyro  dynamics  [5],  Prior  to  introducing 
the  non-linear  system,  Table  1  defines  the  15  state  variables  and  four 
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control  variables.  First,  it  is  important  to  realize  that  the  UAV  will  operate  in  two 
different  sets  of  coordinate  frames:  Body  Frame  (used  for  velocities  and  based 
on  coordinates  with  reference  to  the  actual  vehicle)  and  North-East-Down  Frame 
(used  for  positions  and  based  on  an  initially  fixed  set  of  coordinates  [North  =  x, 
East  =  y,  Down  =  z])  [5], 

Table  1 .  Dynamic  Equation  Variable  Definitions 


Variable 

Variable  Definition 

Unit 

P  x  ’  P  y  1  P  ~ 

Position  Vector  along  the  North-East-Down  (NED)  Frame  (x,y,z) 

meters 

u,v,w 

Velocity  Vector  along  the  Body  Frame  (x,y,z) 

meters/sec. 

<S>,0,¥ 

Roll,  Pitch,  and  Yaw  Angles  (Euler  Angle  in  the  NED  Frame) 

radians 

q,r,s 

Roll,  Pitch,  and  Yaw  Angular  Rates  in  the  Body  Frame 

radians/sec. 

asA 

Longitudinal  and  Lateral  Tip-Path-Plane  (TPP)  Flapping  Angles 

radians 

^ ped,  int 

Intermediate  State  in  Yaw  Rate  Gyro  Dynamics 

N/A 

Pat 

Allows  the  user  to  control  the  ailerons  of  the  Helicopter  UAV 

N/A 

Pan 

Allows  the  user  to  control  the  elevators  of  the  Helicopter  UAV 

N/A 

Pol 

Allows  the  user  to  control  the  collective  pitch  of  the  Helicopter  UAV 

N/A 

5ped 

Allows  the  user  to  control  the  rudder  of  the  Helicopter  UAV 

N/A 

1.  Kinematics 

The  kinematics  component  generates  six  of  the  15  state  variable 
differential  equations  for  our  system.  The  six  equations  are  derived  from  two 
vector  equation  systems,  each  with  three  equations,  based  on  position  and 
velocity  (in  the  x,  y,  and  z  directions).  The  first  equation  set  is  for  position: 

P.=B,-Vb 

where 

Pn=(Px,Py,PZf 
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is  the  position  vector  in  the  NED  Frame  and 


Vb  =  (u,v,w)T 

is  the  velocity  vector  in  the  Body  Frame.  Since  the  position  and  velocity  vectors 
are  defined  in  different  coordinate  spaces,  it  is  imperative  to  have  a 
transformation  matrix  between  the  two.  This  transformation  matrix,  designated 
Bb  ,  is  defined  below  [5], 


bb 


^cos  (9)  cosf't')  sin(<t>)  sin(60  cosfT)  -  cos(<t>)  sinfT)  cos(<t>)  sin(60  cosfT)  +  sin(<t>)  sinfT)  ’ 
cos(60  sinfT)  sin(<t))  sin(60  cosfT)  +  cos(<t>)  cosfT)  cos(<t>)  sin(60  sinfT)  -  sin(<I>)  cosfT) 
-sin(60  sin(<t>)cos(60  cos(4>)cos(60 


The  second  kinematic  equation  set  is  for  rotational  motion: 

=  '^-b 

where 

Q„=(0,£,T)r 

is  the  Euler  Angle  Vector  in  the  NED  Frame  and 

Qh=(q,r,s)T 

is  the  angular  velocity  vector  in  the  Body  Frame.  Just  as  is  the  case  with  the  first 
equation  set,  the  different  coordinate  spaces  require  a  corresponding 
transformation  matrix,  designated  SB  and  defined  below. 


f 


1  tan(#)  sin(cD)  tan(6*)  cos(O) 
0  cos(O)  -sin(O) 

sin(O)  cos(O) 

cos(6*)  cos(6*)  , 


2.  Rigid-Body  Dynamics 

The  six  degrees  of  freedom  of  the  rigid  body  dynamics  of  the  helicopter 
generate  two  additional  systems  of  equations,  with  each  system  containing  three 


9 


equations.  These  six  equations  derive  the  six  components  of  the  helicopter’s 
velocity  (three  in  space  and  three  angular)  and  are  represented  by  the  following 
Newton-Euler  Equations: 


Vb=(-QbxVb)  + 


F  F 

L  (T 


m 


m 


and 


=  r1[M6-(ni)x/.n;,)] 


where 


Fg  =  (mg  sin  6,  mg  sin  (D  cos  6 ,  mg  cos  cD  cos  0)T 
is  the  gravity  force  vector, 


Ft, 

(Xmr  +  XfJ 

Fh  = 

Fby 

= 

A  +  ^  +  y„.  +  n,) 

A. 

(Zmr  +  ZfuS+Zl,f) 

is  the  aerodynamic  force  vector, 


X" 

(A»,-  +  4-/  +Ar) 

Mb  = 

Mby 

= 

1 

Ini 

1 _ 

(Nm.  +  Nvf  +  Nn.) 

is  the  moment  vector,  and 


1  = 


I 

X\ 

0 

0 


0 

r 

y> 

0 


0 

0 
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is  the  moment  of  inertia  matrix.  The  term  m  represents  the  total  mass  of  the 
helicopter  UAV,  and  g  represents  a  constant  value  for  gravity.  The  subscripts 
mr,  tr,  fus,  vf,  and  hf  represent  main  rotor,  tail  rotor,  fuselage,  vertical  fin,  and 
horizontal  fin  or  the  helicopter,  respectively.  The  X,  Y,  Z,  L,  M,  and  N  terms  are 
all  equations  based  on  combinations  of  different  coefficients  and  variables. 
Expanded  detail  of  all  terms  can  be  found  in  [5], 

3.  Main  Rotor  Flapping  Dynamics 


These  two  state  variables,  as  and  bs ,  represent  the  longitudinal  and 

lateral  Tip  Path  Plane  angles  as  the  main  rotor  spins  during  flight.  The  dynamics 
of  these  two  variables  are  represented  by  the  following  equations: 


I 

a.,  =  -r - -  +  — 


f 


da , 


da 


du 


s-u  +  — -w 


dw 


A s 


+  ■ 


Jlon 


d ; 


Ion 


J 


and 


\_5\ 

T  dv 


V  + 


where  r  represents  the  effective  rotor  time  constant  of  the  main  rotor  system, 

das  5bs  das 

T-  and  “7^  are  the  longitudinal  and  lateral  dihedral  effect  derivatives, 

ou  ov  ow 

is  the  flap-back  effect  derivative,  ^Slo„  is  the  effective  linkage  gain  from  the 
control  variable  dlon  to  the  cyclic  pitch  angle  along  the  longitudinal  direction, 

and  ^slal  is  the  effective  linkage  gain  from  the  control  variable  dlat  to  the  cyclic 

pitch  angle  along  the  longitudinal  direction  [5],  Control  variables  are  defined  and 
explained  in  Section  C.  Expanded  detail  of  these  terms  can  be  found  in  [5], 
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4. 


Yaw  Rate  Gyro  Dynamics 


For  small-scale  remote  controlled  helicopters,  the  yaw  moment  is  very 
sensitive;  therefore,  it  is  very  difficult  to  control  yaw  motion  in  manual  flight.  To 
combat  this  challenge,  small-scale  helicopters  are  commonly  equipped  with  a 
yaw  rate  gyro,  which  consists  of  a  gyro  sensor  and  an  embedded  controller,  to 
allow  the  human  controller  to  modify  the  yaw  rate  and/or  heading.  In  modern 
models,  such  a  device  is  not  essential  to  fly  a  programmed  path;  however,  it  is 
reserved  for  the  manual  back-up  flight  system.  Therefore,  the  yaw  rate  gyro’s 
dynamics  must  be  included  in  the  non-linear  flight  model  and  are  represented  by 
the  following  differential  equation: 

^ ped,  int  ^  a  ^ ped  ^ 

where  Ka  represents  a  scaling  value  (a  constant  of  value  3.73),  and  oped 

represents  the  control  variable  for  rudder  servo  input  (defined  and  explained  in 
Section  C)  [5], 

The  15  variables  described  in  the  four  preceding  subsections  are  the 
system  whose  solutions,  based  on  the  variation  of  the  four  control  variables, 
were  used  by  the  team  in  the  non-linear  model  to  code  simulated  flight  paths  for 
the  X-Cell  60  small-scale  UAV  helicopter.  This  same  non-linear  model  is  the 
critical  input  piece  into  the  code  that  will  eventually  create  a  numerical  program  to 
optimize  the  trajectory  of  the  same  UAV  helicopter. 

C.  CONTROLS 

In  order  to  achieve  solutions  to  the  15  state  variable  differential  equations, 
it  is  essential  to  impose  a  set  of  controls  that  will  directly  impact  the  flight  of  the 
helicopter  and,  therefore,  will  generate  state  variable  solutions  throughout  an 
interval  of  time.  These  controls,  which  are  non-dimensional,  have  been  scaled  in 
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this  particular  model  to  have  values  between  -1  and  1.  There  are  four  total 
control  variables,  and  they  are  annotated,  along  with  their  impact  on  the 
helicopter’s  flight,  in  the  following  subsections  [5], 

1.  Aileron  Servo  Input  iplat) 

This  control  variable  allows  the  user  to  control  the  ailerons  of  the 
helicopter.  This  particular  input  will  directly  influence  the  state  of  roll  for  the 
aircraft  (p  and  ®  as  far  as  the  state  variables  are  concerned)  [5], 

2.  Elevator  Servo  Input  (^ton) 

This  control  variable  allows  the  user  to  control  the  elevators  of  the 
helicopter.  This  particular  input  will  directly  influence  the  state  of  pitch  for  the 
aircraft  ( q  and  6  as  far  as  the  state  variables  are  concerned)  [5], 

3.  Rudder  Servo  Input  (dped) 

This  control  variable  allows  the  user  to  control  the  rudder  of  the 
helicopter.  This  particular  input  will  directly  influence  the  state  of  yaw  for  the 
aircraft  ( r  and  ¥  as  far  as  the  state  variables  are  concerned)  [5], 

4.  Collective  Pitch  Servo  Input  (^co/) 

The  collective  pitch  control  is  the  most  unique  of  the  four.  While  it  has  the 
same  range  of  input  values,  it  is  not  an  independent  control.  This  means  that  it 
may  not  be  altered  without  changing  at  least  one  of  the  other  control  variables  to 
balance  it.  Changing  only  the  collective  pitch  will  cause  the  aircraft  to  become 
unstable  and  possibly  crash.  Obviously,  this  provides  disastrous  results  for 
trajectory  simulations  [5], 
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The  user  programmed  inputs  of  the  aforementioned  four  control  variables, 
forced  to  conform  to  pre-determined  constraints,  will  allow  a  property  designed 
non-linear  model  to  create  solutions  to  the  state  variable  differential  equation 
system. 

D.  CONSTRAINT  EQUATIONS 

The  development  of  obstacle  constraints  involves  limitations  to  the 
trajectory  of  the  vehicle  that  do  not  involve  maximum  and  minimum  restrictions 
from  the  system  state  or  control  inputs.  This  section  will  focus  on  obstacle 
avoidance. 

The  main  focus  of  this  particular  UAV  Trajectory  Model  is  for  use  in  an 
urban  environment,  so  it  is  important  to  consider  typical  obstacles  one  might 
encounter  while  travelling  aerially  in  a  city.  There  are  three  obstacle  types: 
buildings,  power  lines  (with  poles),  and  bridges.  Ideally,  a  model  that  can  handle 
all  three  types  is  optimal.  It  is  important  to  note  that  the  model  cannot  produce 
obstacles  with  rigid  corners,  as  this  can  result  in  non-differentiable  functions  as 
we  attempt  to  solve  the  system.  The  obstacle  modeling  in  this  section  uses  the 
p-norm  to  create  a  variation  on  the  standard  3-D  equation  for  an  ellipsoid-shaped 
region.  The  equation  h(x,y,z)  shows  the  general  form  used  to  model  the 
obstacles. 


h(x(t),y(t),z(t))  = 


(  x(t)  -  -  v 


■x„ 


a 


+ 


(  y(t)-ya 


v 


+ 


(  z(t)~ 


-1  =  0 


•  x(t),  y(t),  z(t)  represent  the  position  (x,y,z)  with  respect  to  time 

•  Xp ,  y ,  z represent  the  center  position  locations  inside  of  each  obstacle 

•  a  represents  the  distance  from  the  center  to  the  obstacle  boundary  in 
the  x-direction 

•  b  represents  the  distance  from  the  center  to  the  obstacle  boundary  in 
the  y-direction 
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•  c  represents  the  distance  from  the  center  to  the  obstacle  boundary  in 
the  z-direction  [3], 

The  absolute  value  signs  can  be  removed  from  the  equation  if  the  p- 
values  are  limited  to  even  numbers.  This  is  suitable,  since  increasing  the  value 
of  p  will  merely  alter  the  shape  of  the  obstacle  from  more  circular/elliptical  to 
more  square/rectangular  [3],  Therefore,  the  equation  below,  where  p>0  and 
even,  will  be  used  for  all  subsequent  examples  which  will  model  several 
obstacles  types  using  MATLAB. 


h{x{t),y{t),z{t))  = 


(  x{t)~ 


X„ 


a 


+ 


( y{t)-ym  Y  ( z{t)-z 


+ 


c 


1  =  0 


J 


This  model  assumes  that  the  ground  is  flat  and  represents  the  xy-plane 
which  means  that  the  z-axis  is  vertical  and  z  >  0  in  all  instances. 
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Example  1 
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r  z  —  2 


Y° 
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-1  =  0 


Figure  1 .  A  short,  block  building  with  center  point  (0,0,2)  and  a=4,  b=4,  c=2, 

p=20 
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Example  2 
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15 
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1  =  0 


Figure  2.  A  skyscraper  type  of  building  with  center  point  (0,0,15)  and  a=4, 

b=6,  c=1 5,  p=40 


Example  3 
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20 
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— 1  =  0 
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Bridge  Portion: 
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+ 
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-1  =  0 
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Figure  3.  A  flat  bridge  with  two  end  pillars  (plot  is  constructed  using  three 
separate  functions;  one  for  each  pillar  and  one  for  the  bridge  portion) 


Example  4 
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Suspension  Supports: 


^x  +  15 
v  15-  j 


v°  ( 
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y 

v2y 


+ 


z  + 


2(x  +  15) 


,20 


20 


-1  =  0 


for  x<0 
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,20 


20 
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for  x>0 


Suspension  Cables: 


( x±d v  r r z±(30-c)^ 


20 


+ 


\  .j  j 


+ 


-1  =  0 


\  C  J 


Note:  For  the  Suspension  Cable  Equation,  a  value  of  c  must  be  selected  so  that 
the  z  coordinate  for  the  to  of  the  suspension  cable  matches  the  z  coordinate  for 
the  suspension  support  for  each  selected  value  of  x±d. 


Figure  4.  A  suspension  bridge  with  two  end  pillars  below  the  road  level,  two 
triangular  suspension  supports  above  the  road  level,  and  a  number  of 
cables  (number  and  positions  along  the  x-axis  can  be  altered  easily  by  the 

designer  within  the  equation) 
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Examples  1  through  4  demonstrate  legitimate  urban  obstacle  possibilities 
with  constants  that  the  designer  can  easily  alter  to  change  the  size,  shape,  and 
position  of  the  particular  obstacles.  Although  the  implicit  equation  used  is 
exponential  in  nature,  most  of  the  outcomes  are  fairly  linear  in  nature;  another 
challenge  arises  when  we  concern  obstacles  with  curvature — i.e.,  power  lines. 
This  requires  the  orientation  of  the  previous  model  in  the  direction  of  a  Catenary. 
To  form  the  suspension  supports  in  Example  4,  a  variant  of  the  third  term,  using 
a  linear  x-term  in  addition  to  the  z-term,  provides  the  necessary  result.  Example 
5  demonstrates  an  example  of  a  power  line  obstacle,  but,  contrary  to  Example  4, 
the  x-term  added  to  the  third  term  is  non-linear  in  nature,  which  inherently  results 
in  curvature. 


Example  5 
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1-0  for  c  e[0.5,2] 


Note:  The  larger  the  value  of  c,  the  more  increased  the  curvature  in  the  shape 
becomes.  At  values  of  c>2,  the  plot  in  MATLAB  begins  to  disintegrate.  Values  of 
c>2  should  be  unnecessary  for  the  purpose  of  designing  obstacle  for  this  project. 
For  the  plot  below,  c  =  1.6. 
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Figure  5.  A  set  of  power  lines  with  two  poles  and  a  buffer  area  simulating  the 

curved  area  created  by  hanging  lines 

Examples  1  through  5  have  demonstrated  some  different  equations  for 
modeling  certain  types  of  common  urban  obstacles  into  the  constraint  function. 

E.  THE  PERFORMANCE  MEASURE 

Thus  far,  we  have  discussed  the  development  of  a  non-linear  state 
variable  equation  system,  using  control  variables  and  constraints,  fora  Helicopter 
UAV  model  that  allows  a  user  to  simulate  flights  and  plot  the  trajectories.  A  user 
may  alter  these  trajectories  through  changing  the  control  variables  in  the 
computer  model.  Finally,  an  optimization  function  will  incorporate  the  non-linear 
simulation  model,  and  all  of  its  conditions,  into  an  integral  (evaluated  over  time) 
that  will  allow  the  user  to  minimize  a  particular  aspect  of  the  UAV’s  operation — 
called  a  cost. 

In  order  to  complete  the  problem  formulation  process  for  optimal  control,  it 
is  critical  to  identify  clearly  the  performance  measures  of  the  system  one  wishes 
to  optimize.  A  different  designer  may  choose  different  performance  measures, 
but  for  this  project,  we  will  use  time  and  cost.  An  optimization  function  is  one  that 
will  either  minimize  or  maximize  the  chosen  performance  measures — in  this 
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case,  both  time  and  cost  will  be  minimized.  The  optimization  functions 
developed  in  this  section  will  combine  the  states,  controls,  and  constraints  into 
one  equation  that  eventually  provides  our  desired  optimization  endstate. 

The  problem  in  this  model  finds  a  solution  for  the  optimal  control  vector 
(u*)  which  causes  the  system 


to  follow  an  allowable  trajectory  (x*)  which  minimizes  the  performance  measure 
(J)  such  that 

'/ 

J  =  h(x(tf  ),tf)  + 1  /  (x(/%  u(t),  t)dt 

k 

where  the  function  h(x(tf),tf)  represents  an  endpoint  cost  not  associated  with 
the  trajectory  itself  [6], 

In  this  case,  the  allowable  trajectory  x*  that  minimizes  the  performance 
measure  J  is  called  the  optimal  trajectory  [6],  Two  things  must  be  realized  prior 
to  modeling  optimal  control  problems.  First,  an  optimal  control  solution  may  not 
be  unique.  This  is  not  necessarily  detrimental;  although  it  may  complicate  the 
computational  aspect  of  calculation,  it  can  allow  flexibility  for  the  designer’s 
configuration  scheme.  Second,  an  optimal  control  may  not  exist — this  is 
detrimental  to  our  desired  results 

In  order  for  the  constraint  algorithm  to  perform  properly,  four  critical  input 
components  are  coded  to  formulate  a  solution: 
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•  The  Cost  Function:  Defines  and  codes  the  integral  for  the  actual  cost 
optimization  function 

•  The  Dynamics  Function:  Takes  the  inputted  non-linear  model  from  the 
flight  simulation  codes  to  program  the  state  and  control  differential 
equation  systems  into  the  optimization  algorithm 

•  The  Path  Function:  Allows  the  user  to  program  the  desired  constraints 
into  the  optimization  algorithm  in  order  to  limit  state  variable  outputs 

•  The  Events  Function:  Defines  and  codes  the  initial  and  final  condition 
limits  for  state  variable  values  on  which  the  user  desires  to  enforce 
hard  values. 

The  next  part  of  this  project  will  describe  how  the  algorithm  inside  of  the 
optimization  programs  works  to  produce  a  solution  from  the  aforementioned 
inputs. 

Overall,  this  thesis  will  attempt  to  obtain  a  unique  optimal  solution  for,  at  a 
minimum,  several  examples  of  the  following  scenario  for  the  X-Cell  60  small- 
scale  helicopter  UAV: 

•  Minimize  flight  time;  no  obstacles. 

In  summary,  we  will  find  an  optimal  solution  (minimum  flight  time)  for  J  in 

J  =  j*  /  (x(t),u(t),  t)dt 

o 
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where 
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governed  by 


and 
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III.  PSEUDOSPECTRAL  METHOD  FOR  OPTIMAL  CONTROL 


A.  INTRODUCTION 

For  quite  some  time,  mathematicians  have  struggled  with  a  reliable 
method  for  solving  optimal  control  problems  with  complicated  nonlinear  dynamics 
and  constraints.  Over  the  past  decade,  direct  computational  methods  have 
become  increasingly  popular  for  solving  these  optimal  control  problems  [7],  In 
direct  methods,  a  continuous  time  problem  is  discretized  over  a  set  time  interval, 
and  the  resulting  discrete  problem  is  solved  numerically  using  nonlinear 
programming  algorithms.  Due  to  significant  progress  in  large-scale 
computational  algorithms  and  nonlinear  programming,  direct  Pseudospectral 
(PS)  methods  have  emerged  as  reliable  direct  methods  for  optimal 
control  [8,  9,  10].  For  quite  some  time,  PS  methods  have  been  widely  applied  for 
complex  scientific  computation  models  involving  differential  equations,  and  PS 
methods  have  proven  themselves  as  very  efficient  in  solving  these  problems. 
However,  only  recently  have  PS  methods  intersected  with  nonlinear  optimal 
control  [11,  12].  This  project  uses  a  PS  algorithm  to  numerically  solve  the 
problem  of  optimal  control  for  the  purpose  of  optimal  UAV  trajectory  planning. 

B.  PROBLEM  FORMULATION 

This  section  will  articulate  how  the  PS  method  discretizes  the  problem  of 
optimal  control  subjecting  to  nonlinear  dynamics  into  a  problem  of  finite 
dimensional  nonlinear  optimization,  enabling  the  application  of  computational 
nonlinear  programming  to  solve  for  the  optimal  control  and  optimal  trajectories. 
Consider  the  following  general  optimal  control  problem: 

+i 

min[J(x(»),w(»))]  =  J  F(x(t),u(t))dt  +  E(x(- l),x(l)) 

-i 

subject  to  the  following  set  of  differential  equations  and  initial  conditions: 
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x  =  f  (x,  u)  ■  xer,Mer 


e(x(-l),x(l))  =  0  ;  h{x{t),u{t))  <  0 

where  x  is  the  state  variable  and  u  is  the  control  input  [7], 

E{x{-  l),x(l)) 

is  the  cost  due  to  the  endpoints, 

<?O(-l),x(l))  =  0 

is  the  condition  at  the  endpoints,  such  as  initial  value,  and 

h(x(t),u(t ))  <  0 

is  the  set  of  state  and  control  variable  constraints.  The  PS  method  uses  Sobolev 
spaces  Wm,p ,  which  involve  functions  £(r) ,  defined  in  [— 1, +1] ,  whose  y-th  order 
weak  derivative,  £0) ,  lies  in  the  space  Lp  for  all  0<  j  <m  with  the  norm  [7] 

m 

^  wm’P  ~  | 

7=0 

As  previously  mentioned,  the  PS  method  is  an  efficient  direct  method;  this 
means  that  the  actual  optimal  control  problem,  not  the  associated  necessary 
conditions,  is  discretized  to  obtain  a  non-linear  programming  problem.  Just  as  in 
any  problem  solved  numerically,  the  accuracy  of  the  PS  method  depends 
strongly  on  the  method  of  approximation.  Given  a  function  / (r) : [«, Z?]  — >  IR. ,  one 
could  conventionally  approximate  f(t)  by  interpolating  over  uniformly  spaced 
time  nodes  where  t0  =a,tl  =(b-a)/ N =b,  and  N  is  equal  to  the  number  of 

interpolation  points  (or  time  nodes  in  the  algorithm).  However,  it  is  widely  known 
and  proven  that  uniformly  spaced  points  may  produce  numerical  solutions  with 
much  higher  approximation  error,  for  the  same  number  of  points,  than  those 
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calculated  using  other  more  sophisticated  interpolation  methods  [7], 
Furthermore,  it  is  critical  to  emphasize  that  the  number  of  interpolation  points  in 
calculating  the  solution  to  an  optimal  control  problem  is  not  only  an  issue  of 
efficiency,  but  also  of  feasibility.  A  higher  number  of  interpolation  points  results 
in  a  higher  dimension  in  the  non-linear  programming  model.  If  a  particular  model 
becomes  too  complex,  it  can  easily  overwhelm  computational  capabilities  of  an 
operating  system;  we  obviously  do  not  want  this  to  happen  and  must  carefully 
select  an  efficient  interpolation  method  that  can  provide  a  low-error  solution  with 
relatively  few  nodes  [7], 

The  PS  model  used  for  this  project  involves  interpolating  with  Legendre- 
Gauss-Lobatto  (LGL)  quadrature  nodes.  These  nodes,  denoted  by 
{r0  =  -l<tl<t2<...<tn_l<tn  =1},  are  the  roots  of  LN  where  LN  is  the  A/f/1-order 

derivative  of  the  Legendre  Polynomial  LN(t ) .  Using  this  method,  the  range  of 

integration  is  transformed  universally  to  [-1,+1],  which  is  the  interval  for  Legendre 
Polynomials.  Although  the  LGL  interpolation  points  are  not  evenly  spaced,  they 
are  symmetric  about  the  midpoint  0  [7],  Figure  6  shows  a  range  of  LGL  points  for 
N  =  16. 


i - 1 - 

LGL  Points;  N=16 
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Figure  6.  LGL  Nodes,  N  =  16 
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PS  methods  are  widely  applied  to  numerically  solve  models  using  partial 
differential  equations  (PDEs).  However,  optimal  control  problems  have  several 
fundamental  differences  from  the  computation  of  PDEs.  Solving  optimal  control 
problems  asks  for  the  collective  and  simultaneous  solving  of  several  different 
systems,  including  the  differential  equation  governing  the  control  system,  the 
integration  of  the  cost  function,  and  the  state-control  constraints.  Then  these 
approximations  are  integrated  together  to  form  a  problem  of  discrete  nonlinear 
optimization  which  must  be  solved  numerically  to  find  the  approximate  optimal 
control  [7], 

In  a  PS  optimal  control  method,  the  state  and  control  functions,  x(t)  and 
u(t),  are  approximated  by  Nth  order  Polynomials  based  on  interpolation  at 
selected  LGL  quadrature  nodes.  In  the  discretization,  the  state  variables  are 
approximated  by  the  vectors 


xNk 


v 


j 


which  is  an  approximation  of  x(tk) .  Also,  the  vector  uNk  is  the  approximation  of 
u(tk) .  Consequently,  a  discrete  estimation  of  the  function  xt(t)  is  the  vector 

—  N 1  —N  2  —NN~\ 

xf  ,  X.  J. 


A  continuous  estimate  can  be  defined  by  a  polynomial  interpolation  xtN(t)  where 


k= 0 


where  0A(r)  is  the  Lagrange  interpolating  polynomial  [13].  In  contrast,  the 
control  input  is  approximated  by  the  non-polynomial  interpolation 
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xrN(t)-f(xN(t)) 
g(xN(t )) 


In  the  notations,  the  discrete  variables  are  denoted  by  letters  with  an  upper  bar, 
such  as  xtm  and  uNk .  If  k  in  the  superscript  and/or  /  in  the  subscript  are  missing, 

then  the  notation  represents  a  vector  or  matrix  in  which  those  particular  indices 
run  from  their  minimum  to  maximum  values  [7],  For  example, 


—  Wl  —N  2 
Xi  ,  X.  , 


-NN 


X 


Nk 


v 
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( - 

-WO 

—  Wl 

—  NN  \ 

X\ 
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■"  ^1 
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—  Wl 
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x9 
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—  Wl 

—  NN 

—  \ 
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— W  2 
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,u  ,  • 

~>U  J 

where 


is  the  interpolation  point  (node)  number 


i  —  { js  the  state  variable  number 

N  =  the  total  number  of  interpolation  points  (nodes) 
r  =  the  total  number  of  state  variables 
and  k,  i,  N,  and  rare  the  subscript  and  superscript  definitions  [7], 
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The  analysis  of  spectral  methods  demonstrates  that  PS  method  is  simple, 
accurate,  and  relatively  fast  in  the  estimation  of  smooth  functions,  integrations, 
and  differentiations.  These  are  all  critical  pieces  for  solving  optimal  control 
problems.  The  derivative  of  xtN{t)  at  the  LGL  node  tk  is  easily  computed  by 
using  matrix  multiplication: 


*"(*0 )*/V(0*/V(?2)-"*/V(^)  *  =  D(XD 


N  . 


.  N 


N  , 


<—N\T 


where  the  (jV  +  l)x(JV  +  l)  differentiation  matrix  D  is  defined  by 


A*  =< 


( Ml 

v  (tic ) )  V  t  ~  4  j 


,if  i^k 


N(N  +  \ )  ...  ,  . , 

-,if  i-k-N 


4 

0, 


otherwise 


The  cost  functional  J[x(-),w(-)]  is  approximated  by  the  Gauss-Lobatto  integration 
rule, 

J [*(•) , «(•)] * 7"  (F ¥ , u" )  =  £ F (xm , u‘ * )Wl+E(xm, xm ) 

A-=0 


where  wk  are  the  LGL  quadrature  weights  defined  by 

-  1  f  2  1 

W*  =  (A(4))2  vJV(JV+i)J’ 

The  approximation  is  so  accurate  that  it  has  no  error  if  the  integrand  function  is  a 
polynomial  of  a  degree  less  than  2N  [7], 
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In  order  to  demonstrate  the  PS  discretization,  consider  the  following 
general  nonlinear  optimization  problem: 

Find  xNk  ef  and  uNk  eM,  k  =  {0,1,...,  TV) ,  that  minimize 
7N(xN,uN)=YJF(xNl,  uNl  K  +  E(xm,xNN ) 

A'=0 


subject  to 


Zn  — Ni  /•/ — Nk  — Nk  \  e 

Dkix  -f(x  ,u  )  <8 

i=0  oo 

e(xN\xm )  <8 

00 


h(xm ,um)<  S 

where  8  is  a  small  positive  number  as  the  relaxation  of  the  constraints.  Such 
relaxation  is  necessary  as  shown  in  [14],  It  guarantees  the  feasibility  of  the 
problem. 

N 

Zn  — Ni  r( — Nk  — Nk\  ^  c 

Dklx  ~f(x  ,u  )  <8 

i— 0  oo 

represents  the  discretization  of  the  control  system  defined  by  the  differential 
equation.  Constraints  are  imposed  and  then  tightened,  as  necessary,  to  define  a 
search  region  within  which  nonlinear  programming  software  packages  can 
produce  a  feasible  optimal  control  solution  [7],  In  [7],  it  is  proven  for  feedback 
linearizable  systems — a  family  of  widely  used  control  systems — that  the  value  of 
the  optimal  cost  of  the  discretized  nonlinear  optimization  converges  to  the 
optimal  cost  of  the  original  problem  of  optimal  control.  It  is  also  proven  that  the 
method  has  a  high  order  rate  of  convergence  depending  on  the  smoothness  of 
the  original  problem.  This  result  implies  that  the  PS  methods  are  convergent  and 
numerically  efficient. 
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In  conclusion,  this  chapter  has  shown  how  PS  optimal  control  methods 
allow  a  problem  of  optimal  control  subject  to  complex  nonlinear  dynamical  and 
algebraic  constraints  to  be  discretized  and  solved  numerically.  The  discretization 
works  in  a  harmonic  way  for  the  multiple  components  in  the  problem,  the  cost 
function,  the  nonlinear  control  system,  and  the  algebraic  constraints.  The 
efficient  transition  from  continuous  to  discrete  is  crucial  in  the  solution  of  the  UAV 
optimal  trajectory  problem.  Also,  further  analysis  of  References  7,  14  and  15 
show  that  calculating  optimal  cost  solutions  using  PS  methods  has  a  very  high 
rate  of  convergence  when  compared  to  other  methods;  in  fact,  the  Legendre  PS 
method  will  have  a  faster  convergence  rate  than  any  polynomial  method  [7,  14, 
15]. 
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IV.  SIMULATIONS 


A.  INTRODUCTION 

The  objective  for  this  thesis  is  to  computationally  find  optimal  trajectories 
and  controls  and  obtain  a  unique  solution  for,  at  a  minimum,  the  scenario 
discussed  at  the  end  of  Chapter  II:  minimize  flight  time  with  no  obstacles. 

Before  simulating  this  scenario,  it  is  necessary  to  mesh  two  existing  codes 
from  two  separate  entities  in  order  to  use  an  optimizing  algorithm  in  conjunction 
with  the  non-linear  model  of  the  X-Cell  60  Helicopter  UAV.  Prior  to  meshing 
these  codes  together,  it  is  critical  to  ensure  that  they  both  function  separately. 
Consequently,  flight  simulations  are  conducted  using  the  non-linear  helicopter 
and  flight  simulation  codes  provided  by  the  team  from  the  National  University  of 
Singapore,  and  a  simple  example  scenario  is  optimized  using  an  algorithm 
developed  at  the  Naval  Postgraduate  School.  Once  the  codes  operate 
separately,  the  non-linear  helicopter  model  is  imported  into  the  optimizing 
program  to  provide  a  single  code  that  will  optimize  the  trajectory  for  our  selected 
Helicopter  UAV  model. 

B.  FLIGHT  SIMULATION 

The  codes  used  from  the  National  University  of  Singapore  contain  two  key 
files:  1)  The  nonlinear  helicopter  model  and  2)  The  flight  simulation  code  (which 
takes  the  non-linear  model  and  inputs  it  into  a  fourth  order  Runge-Kutta  Method 
loop  to  fly  over  a  set  number  of  time  steps).  Initially,  the  code  designers  set  the 
control  variables  8lat ,  Slon ,  and  Scol  as  fixed  values  throughout  all  time  steps  and 

the  8  d  control  variable  to  change  linearly  with  respect  to  the  Intermediate  State 

in  Yaw  Rate  Gyro  Dynamics  state  variable.  These  control  variable  settings 
ensured  flight  stability  and  are  set  in  a  manner  to  fly  the  UAV  on  a  relatively 
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straight  path  in  the  y  direction  (due  East).  Figures  7  and  8  show  the  positions 
and  velocities,  respectively,  in  the  x,  y,  and  z  directions  over  time  with  the 
parameter  set  by  the  code  designers. 


Position  described  in  ground  frame 


Figure  7.  Flight  Simulation,  positions  (x,y,z)  with  respect  to  time 


Velocity  described  in  body  frame 


Figure  8.  Flight  Simulation,  velocities  (u,v,w)  with  respect  to  time 
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As  shown  in  Figures  7  and  8,  the  Flight  Simulation  code,  using  the  non¬ 
linear  model,  produces  a  stable  flight  for  the  Helicopter  UAV;  therefore,  the  non¬ 
linear  model  should  behave  properly  when  integrated  into  the  optimization 
algorithm. 

C.  OPTIMIZATION  EXAMPLE 


In  order  to  gain  a  better  understanding  for  the  optimization  program,  it  was 
essential  to  run  an  elementary  example  prior  to  solving  such  a  complex  system 
as  our  Helicopter  UAV  model.  The  initial  example  is  a  simple  system  with  two 
state  variables  and  one  control  variable.  The  example  problem  is 

i  =  T 

< 

y  =  y  +  u 

where  u  is  the  control  variable  and  the  performance  measure  is 


1 

min  [  J]  =  |  (x2  +  y2  +  .005 u)  dt 

0 


subject  to  the  constraint 


and  the  initial  conditions 


T-8 


f  i  V 


t  — 

v  2  j 


2 


x(0)=0 
^(0)  =  -l 


where  t  e  [0,1] . 


It  turns  out  that  the  optimization  algorithm  provides  an  optimal  solution  to  this 
problem.  This  is  determined  by  the  Exit  Codes;  these  codes  are  built  into  the 
program  to  show  that  a  solution  is  optimal  or  to  help  refine  possible  reasons  for 
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non-optimal  solutions.  In  this  particular  case,  the  Exit  Code  is  1 — an  optimal 
solution  with  all  conditions  satisfied  [16].  Figure  9  depicts  the  control  variable 
values  throughout  the  time  interval,  and  Figure  10  shows  the  x,  y,  and  constraint 
values  throughout  the  time  interval.  It  is  clear  from  Figure  10  that  the  state 
variable  values  follow  the  constraint.  All  of  the  outputs  from  the  optimization 
algorithm  determine  that  the  value  J  =  .1803  is  a  minimum  value. 


Plot  of  control  variable  value  at  time  (t) 


Figure  9.  Control  Variable  (u)  values  at  (t)  for  the  Optimization  Example 
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Plot  of  x,y,and  constraint  values  at  time  (t) 


Figure  10.  x,  y,  and  constraint  values  at  (t)  for  the  Optimization  Example 

D.  OPTIMAL  CONTROL  OF  THE  HELICOPTER  UAV 

Since  both  codes  have  separately  provided  successful  results,  it  is  time  to 
mesh  the  codes  and  begin  working  through  some  examples  using  our  non-linear 
Helicopter  UAV  model.  It  is  crucial  for  us  to  understand  that  simulations  are 
testing  the  optimal  control  algorithm  and  the  non-linear  UAV  model.  The  UAV 
model  has  been  verified  for  stability  and  simulation  of  smooth  curve  flight  but  not 
necessarily  for  use  in  optimal  control.  Any  simulation  errors  or  inconsistencies 
could  result  from  incompatibilities  in  either  code. 

For  all  optimization  simulations,  30  LGL  time  nodes  are  used.  This 
number  of  nodes  seem  to  provide  the  best  results  throughout  simulations.  The 
control  variables,  as  given  in  the  non-linear  helicopter  model,  have  non- 
dimensional  values  and  are  bounded  in  the  following  manner: 
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3.,  6  [-U] 

4^[-i.i] 

8Co,^l- 1.0] 

<^6[-l,l]. 

1.  Minimize  Time,  No  Obstacles 

For  the  time  minimization  performance  measure,  four  simple  examples  are 

run: 

Example  1:  The  UAV  travels  from  the  point  (0,  0,  -10)  to  the  point  (10,  0,  -10) 
(the  initial  and  final  position  values  in  the  y  and  z  directions, 
respectively,  are  the  same) 

Example  2:  The  UAV  travels  from  the  point  (0,  0,  -10)  to  the  point  (10,  10,  -10) 
(the  initial  and  final  position  values  are  the  same  in  only  the  z 
directions) 

Example  3:  The  UAV  travels  from  the  point  (0,  0,  -10)  to  the  point  (50,  0,  -10) 
(same  conditions  in  the  y  and  z  directions  with  the  distance 
travelled  in  the  x  direction  increased  from  10  to  50  units) 

Example  4:  The  UAV  travels  from  the  point  (0,  0,  -10)  to  the  point  (50,  50,  -10) 
(the  total  distance  travelled  in  both  the  x  and  y  directions  is 
increased  from  10  to  50  units). 

All  four  examples  have  the  same  initial  conditions  for  the  state  variables  (all 

values  are  zero,  except  for  pz ,  which  is  initially  -10;  this  implies  that  the 

helicopter  begins  at  a  position  10  units  above  the  ground).  Since  the  objective  is 

to  minimize  time,  the  integral  for  the  performance  measure  is  simply 
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min  [  J]  =  j*  1 

0 

where  the  solution  simplifies  to  tf . 

For  Example  1,12  of  the  state  variables  have  fixed  boundary  conditions  at 
tf ;  they  are  annotated  in  the  vector  equation 

10  " 
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-10 
0 
0 
0 
0 

o  ■ 

0 
0 
0 

0  , 


Since  there  are  no  obstacles  impeding  the  trajectory,  the  expectation  is 
that  the  solution  should  be  a  relatively  straight  flight  path  in  the  x-direction  with  a 
smooth  acceleration  for  the  first  half  of  the  flight  followed  by  a  smooth 
deceleration  for  the  second  half  in  order  to  return  to  a  u  and  v  value  of  0  at  tf . 

As  shown  in  Figures  11  and  12,  the  actual  solution  provided  by  the  algorithm 
adheres  closely  to  our  predictions. 
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Positions  over  Time  for  Example  1  (min.  time,  no  constraints;  x(0)  =  0,  x(tf)  =  10) 


Figure  1 1 .  Positions  over  time  for  Example  1 


Velocities  over  Time  for  Example  1  (min.  time,  no  constraints;  x(0)  =  0,  x(tf)  =  10) 


Figure  12.  Velocities  over  time  for  Example  1 


The  minimum  time  required  is 


min  [J]  =  tf  =2.6905  . 
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Figure  13  depicts  the  four  control  variable  values  throughout  the  flight  time 
span.  In  the  control  variable  plots  for  all  examples  in  this  chapter,  the  following 
substitutions  apply: 


Ul-S, 


lat 


U2  =  S, 


Ion 


m=s 


col 


u*=spei. 


Control  Variable  Values  over  Time  for  Example  1  (min.  time,  no  constraints;  x(0)  =  0,  x(tf)  =  10) 


Figure  13.  Control  variable  values  over  time  for  Example  1 

The  solution  value  and  the  position  plots  in  Figure  1 1  seem  to  hint  that  our 
solution  is  valid  and  quite  possibly  a  true  optimal  solution  to  Example  1. 
However,  the  velocity  and  control  variable  plots  show  some  unexpected 
inconsistencies  that  contradict  the  Exit  Code  from  the  program  output  (which 
suggests  that  the  solution  is  truly  optimal).  The  Exit  Code,  1,  as  defined  in  the 
algorithm’s  user  guide,  tells  us  that  the  simulation  has  produced  an  optimal  value 
for  the  performance  measure  [16],  Figure  12  shows  an  unexpected  range  in 
velocities  in  the  z  direction — there  should  be  very  little  velocity  in  the  y  or  z 
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direction  to  fly  from  0  to  10  in  only  the  x  direction.  Also,  Figure  13  demonstrates 
high  levels  of  instability  in  the  control  variable  values  near  t0  and  tf.  These 

observations  from  the  simulation  results  lead  us  to  believe  that  the  trajectory 
does  not  obey  the  dynamical  rules  of  the  system.  This  issue  will  be  discussed 
later  in  this  section. 

For  Example  2,  12  of  the  state  variables  have  fixed  boundary  conditions  at 
tf ;  they  are  annotated  in  the  vector  equation 
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Since,  again,  there  are  no  obstacles  involved,  the  expectation  is  that  the 
solution  should  be  a  relatively  straight  flight  path  in  the  xy-direction  with  a  smooth 
acceleration  for  the  first  half  of  the  flight  followed  by  a  smooth  deceleration  for  the 
second  half  in  order  to  return  to  a  u  and  v  value  of  0  at  tf .  One  would  expect  a 

similar  position  and  velocity  curve  in  the  x  and  y  directions.  Figures  14  and  15 
show  that  the  actual  solution  provided  by  the  algorithm  adheres  closely  to  our 
predictions. 


44 


Positions  over  Time  for  Example  2  (min.  time,  no  constraints;  x(0),y(0)  =  (0,0),  x(tf),y(tf)  =  (10,10) 


Figure  14.  Positions  over  time  for  Example  2 


Velocities  over  Time  for  Example  2  (min.  time,  no  constraints;  x(0),y(0)  =  (0,0),  x(tf),y(tf )  =  (10,10) 


Figure  1 5.  Velocities  over  time  for  Example  2 
The  minimum  time  required  is 


min  [J]  =  tf  =2.7255  . 
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Figure  16  depicts  the  four  control  variable  values  throughout  the  flight  time 


span. 


Control  Variable  Values  over  Time  for  Example  2  (min.  time,  no  constraints;  x(0),y(0)  =  (0,0),  x(tf),y(tf)  =  (10,10) 


Figure  16.  Control  variable  values  over  time  for  Example  2 

The  solution  value,  the  position  plot  in  Figure  14,  and  the  velocity  plot  in 
Figure  15  (in  the  x  and  y  directions)  all  seem  to  hint  that  our  solution  is  valid  and 
quite  possibly  a  true  optimal  solution  to  Example  2.  However,  the  velocity  in  the 
z  direction  and  control  variable  plots  show  some  unexpected  inconsistencies  that 
help  to  explain  the  Exit  Code  from  the  program  output  (which  suggests  that  the 
solution  is  not  optimal).  Example  2  had  a  different  Exit  Code  output  from 
Example  1  (Exit  Code  41)  which,  according  to  the  user’s  guide,  means  that  our 
solution  is  nearly  optimal  but  encounters  numerical  difficulties  when  the  sub- 
optimal  solution  cannot  be  improved  [16].  This  is  apparent  from  the  simulation 
run  times:  10  to  15  seconds  for  Exit  Code  1  but  400  to  1300  seconds  for  Exit 
Code  41.  Similar  to  Example  1,  Figure  15  shows  an  unexpected  range  in 
velocities  in  the  z  direction — there  should  be  very  little  velocity  in  the  z  direction 
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to  fly  from  0  to  10  in  the  x  and  y  directions — and  Figure  16  demonstrates  high 
levels  of  instability  in  the  control  variable  values  near  t0  and  tf  which  could 

directly  influence  the  optimality  of  the  solution. 

For  Example  3,  12  of  the  state  variables  have  fixed  boundary  conditions  at 
tf ;  they  are  annotated  in  the  vector  equation 


Similar  to  Example  1 ,  the  UAV  should  fly  a  relatively  straight  path  in  the  x- 
direction  with  a  smooth  acceleration  for  the  first  half  of  the  flight  followed  by  a 
smooth  deceleration  for  the  second  half  in  order  to  return  to  a  u  and  v  value  of  0 
at  tf.  As  shown  in  Figure  17,  the  actual  solution  provided  by  the  algorithm 

adheres  closely  to  our  predictions.  However,  Figure  18  shows  unexpected 
changes  in  velocity  in  both  the  y  and  z  directions. 
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Position  over  Time  for  Example  3  (min.  time,  no  constraints;  x(0)  =  0,  x(tf)  =  50) 


Figure  17.  Positions  over  time  for  Example  3 


Velocities  over  Time  for  Example  3  (min.  time,  no  constraints;  x(0)  =  0,  x(tf)  =  50) 


Figure  18.  Velocities  over  time  for  Example  3 

Figures  17  and  18  clearly  do  not  match  up  for  their  y  and  z  components. 
Based  on  Figure  18,  there  should  be  substantial  position  changes  in  all  three 
directions  in  Figure  17.  This  is  not  the  case;  however,  the  simulation  achieves  an 
Exit  Code  1  value  that  implies  the  following  is  an  optimal  solution  [16]: 
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min  [J\=tf  =5.231 1 . 

However,  the  plot  discrepancies  suggest  an  issue  with  the  non-linear  UAV  model 
portion  of  the  code.  Figure  1 9  demonstrates  similar  endpoint  instability  for  the 
control  variable  as  were  seen  in  Examples  1  and  2;  however,  the  instability  is  not 
nearly  as  drastic  as  in  the  first  two  examples. 


Control  Variable  Values  over  Time  for  Example  3  (min.  time,  no  constraints;  x(0)  =  0,  x(tf)  =  50) 


Time  (t) 

Figure  19.  Control  variable  values  over  time  for  Example  3 

For  Example  4,12  of  the  state  variables  have  fixed  boundary  conditions  at 
tf ;  they  are  annotated  in  the  vector  equation 
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Since,  again,  there  are  no  obstacles  involved,  the  expectation  is  that  the 
solution  should  be  a  relatively  straight  flight  path  in  the  xy-direction  with  a  smooth 
acceleration  for  the  first  half  of  the  flight  followed  by  a  smooth  deceleration  for  the 
second  half  in  order  to  return  to  a  u  and  v  value  of  0  at  tf .  Just  as  in  Example  2, 

one  would  expect  a  similar  position  and  velocity  curve  in  the  x  and  y  directions. 
Figures  20  and  21  show  that  the  actual  solution  provided  by  the  algorithm 
adheres  closely  to  our  predictions. 
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Positions  over  Time  for  Example  4  (min.  time,  no  constraints;  x(0),y(0)  =  (0,0),  x(tf),y(tf)  =  (50,50) 


Figure  20.  Positions  over  time  for  Example  4 


Velocities  over  Time  for  Example  4  (min.  time,  no  constraints;  x(0),y(0)  =  (0,0),  x(tf),y(tf)  =  (50,50) 


Figure  21 .  Velocities  over  time  for  Example  4 
The  minimum  time  required  is 


=6.3853 . 
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Figure  22  depicts  the  four  control  variable  values  throughout  the  flight  time 


span. 


Control  Variable  Values  over  Time  for  Example  4  (min.  time,  no  constraints;  x(0),y(0)  =  (0,0),  x(tf),y(tf)  =  (50,50) 


Figure  22.  Control  variable  values  over  time  for  Example  4 

Just  as  in  Example  2,  the  solution  value,  the  position  plot  in  Figure  20,  and 
the  velocity  plot  in  Figure  21  (in  the  x  direction)  all  seem  to  hint  that  our  solution 
is  valid  and  quite  possibly  a  true  optimal  solution  to  Example  4.  However,  the 
velocity  in  the  y  and  z  directions  and  control  variable  plots  show  some 
unexpected  inconsistencies  that  help  to  explain  the  Exit  Code  of  41  from  the 
program  output  (which,  just  as  in  Example  2,  suggests  that  the  solution  is  not 
optimal)  [16].  The  inconsistencies  between  the  position  and  velocity  plots  are 
very  similar  to  those  in  Example  3.  The  control  variable  values  in  Figure  22  show 
that  the  endpoint  instability  is  not  nearly  as  drastic  as  in  Examples  1  and  2  but 
similar  to  that  shown  in  Example  3.  It  seems  that  lengthening  that  path  tends  to 
smooth  the  control  variable  variation  to  some  degree. 

The  previous  four  simple  examples  have  demonstrated  that  the  meshing 
of  the  nonlinear  UAV  model  and  the  optimization  algorithm  codes  is  functional 
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and  effectively  simulates  a  three  dimensional  trajectory  while  outputting  an 
optimal  or  near  optimal  solution  for  minimal  flight  time  for  several  different  sets  of 
initial  and  final  conditions.  However,  there  are  several  issues  that  must  be 
worked  out  in  the  future.  The  inconsistency  between  the  path  and  velocity 
trajectories  implies  that  discrepancies  exist  when  the  UAV  model  is  integrated 
with  the  optimal  control  software  package.  It  is  a  known  fact  that  computational 
optimal  control  requires  higher  accuracy  than  what  is  needed  in  the  integration  of 
ODEs.  As  recently  as  two  weeks  ago,  the  UAV  Laboratory  of  the  National 
University  of  Singapore  developed  a  new  model  of  the  same  UAV  system. 
According  to  their  analysis  and  simulations,  the  accuracy,  when  comparing  to 
experimentation  data,  is  much  higher  then  the  old  model  that  is  being  used  in  our 
program.  We  are  advised  to  use  the  new  model  in  all  simulations.  The  next  step 
of  research  is  to  integrate  the  new  model  into  the  package  developed  in  this 
thesis  to  achieve  improved  simulation  results.  In  addition,  various  numbers  of 
nodes  and  scenarios  are  to  be  numerically  tested  and  analyzed. 
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V.  CONCLUSIONS  AND  FUTURE  RESEARCH 


A.  CONCLUSIONS 

This  thesis  project  has  been  a  very  interesting  look  into  the  value  of 
numerical  methods  for  solving  complex  problems.  There  are  several  findings, 
some  positive  and  some  negative,  that  have  impacted  the  results  and, 
consequently,  have  set  the  ground  work  for  a  plethora  of  future  research. 

First,  and  foremost,  the  meshing  of  two  models  from  two  different  software 
packages  into  one  set  of  optimization  codes  was  effective  and  produced  optimal 
results  for  several  simple  examples.  The  organization  of  the  optimization  codes 
into  five  separate  user  input  file — one  for  the  performance  measure,  one  to 
manage  the  system  of  dynamic  and  control  equations,  one  to  manage  the 
constraints,  one  to  manage  the  initial  and  final  conditions  of  the  equation  system, 
and  one  main  code  to  run  the  algorithm — proved  very  useful  in  segregating  and 
resolving  errors  during  the  early  coding  phases.  Also,  this  organization  of  the 
codes  made  it  quite  simple  to  set  up  new  examples  for  analysis. 

Interestingly,  the  optimization  code  is  very  sensitive  to  bounds  placed  on 
the  variables.  It  was  an  obvious  assumption  that  bounds  too  restrictive  in  nature 
would  make  the  problem  infeasible.  As  a  result,  the  program  cannot  calculate  an 
optimal  solution.  However,  if  the  range  between  the  upper  and  lower  bound  is 
too  large,  it  may  take  a  long  time  for  the  program  to  converge.  In  worse  cases, 
the  program  may  diverge.  This  sensitivity  forces  the  user  to  carefully  select  the 
bounds  and  continue  refining  them,  in  some  cases,  to  allow  the  program  to  run 
successfully. 

Since  some  examples  resulted  in  an  optimal  solution  and  some  did  not,  it 
seems  logical  that  there  are  some  issues  with  both  models  that  need  to  be 
worked  out  in  order  to  make  the  optimization  algorithm  and  nonlinear  model 
operate  in  harmony.  Since  the  non-linear  model  is  not  originally  designed  for 
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optimization,  this  is  a  sensible  conclusion.  Some  of  the  constant  values  that 
create  an  ideal  environment  for  flight  simulation  may  not  necessarily  do  so  for 
optimization  of  a  flight  trajectory.  The  fact  that  it  works  in  certain  cases  is  a  sign 
that  the  nonlinear  model  should  not  have  to  be  altered  much  to  operate  much 
more  effectively  within  the  optimization  algorithm  for  additional  simulations. 

Another  disparity  that  hints  toward  non-conformity  within  the  nonlinear 
UAV  model  is  the  issue  of  disagreement  between  one  or  more  of  the  velocity 
components  over  time  as  compared  to  their  respective  position  components. 
The  fact  that,  in  some  cases,  the  optimization  algorithm  outputs  an  Exit  Code  1 
(an  optimal  solution  output)  while  the  plots  do  not  make  sense,  hints  that  there  is 
a  discrepancy  when  the  UAV  model  is  integrated  with  the  optimal  control 
program.  These  inconsistencies  must  be  adjusted  before  one  can  simulate  and 
optimize  more  complex  trajectories. 

The  final  finding,  which  would  make  the  examples  in  this  thesis  very 
difficult  to  transition  to  real-life  flight  simulations,  is  the  erratic  nature  of  the 
control  variable  values,  especially  near  the  initial  and  final  time  steps.  It  seems 
that  the  actual  helicopter  UAV  would  fly  in  a  very  unstable  fashion  with  the 
plotted  control  variable  inputs.  Smoothing  these  input  values  would  make  the 
model  a  more  realistic  feat  for  actual  test  flight  trajectory  analysis. 

Overall,  this  project  was  an  interesting  and  enjoyable  research  experience 
where  I  was  able  to  learn  a  great  deal  about  a  subject  of  which  I  previously  knew 
almost  nothing  about.  The  meshing  of  two  codes  into  one  functional  optimization 
algorithm  with  the  ability  to  simulate  an  extremely  complex  flight  trajectory  and 
output  an  optimal  solution  was  a  great  success,  and  I  look  forward  to  pursuing 
future  endeavors  on  this  project  with  my  advisor. 
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B.  FUTURE  RESEARCH 


There  is  a  tremendous  amount  work  that  can  be  done  to  further  the 
research  in  this  thesis,  and  I  plan  on  pursuing  some  of  this  work  with  Dr.  Kang  in 
the  future.  Some  future  endeavors  include  the  following: 

•  Working  with  the  recently  improved  UAV  model  to  get  it  to  conform 
better  to  the  optimization  algorithm 

•  Applying  the  optimization  algorithm  with  this  particular  helicopter 
UAV  to  more  complex  flight  trajectories  such  as  collision  avoidance 
in  obstacle  dense  environment 


•  Optimizing  other  performance  measures 

•  Applying  the  same  optimization  algorithm  to  other  UAV  types  and 
actual  helicopters. 


The  modeling  of  these  scenarios  could  prove  invaluable  to  the  United  States 
Military’s  ability  to  maximize  the  effectiveness  of  limited  UAV  resources  in 
combat  zones  and  could  eventually  model  trajectories  minimizing  the  risk  to 
pilots  and  aircraft  for  helicopters  of  different  size  and  purpose. 


Improving  both  the  UAV  model  and  the  optimization  algorithm  will  enable 
the  program  to  calculate  solutions  for  more  complex  trajectory  simulations.  For 
example,  a  user  could  add  obstacles  to  the  UAV’s  path,  which,  up  to  this  point, 
we  have  not  been  able  to  code  effectively.  Also,  one  could  track  a  path 
(potentially  a  hostile  person  or  vehicle  on  the  ground)  with  minimum  variance 
from  that  path.  This  model  would  be  incredibly  useful  in  urban  environments.  In 
addition  to  minimizing,  a  user  could  opt  to  try  to  maximize  a  UAV’s  capability, 
such  as  the  number  targets  or  amount  of  ground  area  observed  within  a  given 
flight  time  window.  In  the  examples  for  this  thesis,  we  minimized  the  flight  time, 
but  additional  parameters  could  be  added  to  the  model  to  minimize  fuel 
consumption  or  flight  cost.  These  are  just  a  few  examples  of  more  complex 
trajectories  with  alternate  performance  measures  for  realistic  UAV  operations. 
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Trajectory  optimization  for  helicopter  UAVs  definitely  has  real-world 
applications,  which  means  many  more  ways  to  continue  building  on  this  project. 
Hopefully,  these  applications  will  provide  the  military  with  plausible  options  for 
maximizing  effectiveness  and  minimizing  risk  and  cost  on  the  future  battlefield. 
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